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Abstract 

The problem of determining 3D density fields from single 2D projections is 
hopelessly under-determined without additional assumptions. While parame- 
terized inversions are typically used to solve this problem, we present theoretical 
results along a different route to the elimination of indeterminacy. More specif- 
ically, we consider the case in which radiography is being used to study objects 
or processes evolving over some time interval. This evolution can be measured 
at several points, generating a sequence of radiographs. A priori knowledge of 
constraints on possible evolution permits us to rule out sequences which are not 
consistent with those constraints. If we now consider the space of possible mea- 
surement sequences to be the data space, we find that, under the influence of 
the constraints, our data space has increased in dimension while the dimension 
of the space of unknowns has remained the same. When enough measurements 
have been made, inversion of a dynamically constrained sequence of single angle 
radiographs becomes possible. 
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1 Introduction 



We outline the paper and briefly introduce the ideas motivating the results we obtained 



2 The Problem 



Our particular problem, that of deducing 3-dim density fields from 2-dim radiographic pro- 


jections of that field, is introduced here. The abstract formulation is then given. We arrive 


at a mathematical problem arising in many 


data driven analysis situations. Predictably, it 


goes by different names in different contexts. 





3 Notation 



Oh, the pleasures and miseries of the use and abuse of creative naming and labeling! Wej 



try to generate a notation that helps the readers in their quest to follow the derivations. We| 
resist the temptation to make our notation look good by pointing out examples of notoriousj 
notation 6 



4 Linear Stationarity implies easy solutions 



That is, conceptually easy ... easy might be very difficult indeed if the dimensions involved 



are big enough! In this section we reduce our problem to controlling the dimension of a 



certain sequence of intersections (which we want to force to 0!) 



5 Trans versality is (more than) enough 



What conditions on those intersections, which we met in the previous section, guarantee that 



the dimension of those intersections goes to 0? Transversality does! In fact it implies that 



the convergence to is as fast as possible! We include a result about how slow the dimension 



can go to zero and still get there 



6 Optimal Dynamics are generic 



Does a typical problem lead to a sequence of intersections whose dimensions converge to zero 



as fast as possible? The answer is that in a carefully defined sense YES 
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7 Extension to the Nonlinear Case 



So far everything we have developed has been for the case in which the underlying dynamics] 



are linear. What can we say about the nonlinear case? We use our hands (waving then 
about) to say something which we make respectable by formulating two conjectures anc 



sketching how we think it should go. 19 
8 Relation to Known Results] 



Of course even though this work was new to us (and therefore fun!) there is nothing nc 



under the sun and we find other work that is very similar. Our hope is that the difference; 



will be useful to the readers and that the readership will be much wider than the few experts] 



who would find the contents not surprising 



20 



9 Numerical Examples 



Here we demonstrate the observability of a particular linear system as well as some of the 
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busy most researchers are and how they tend to read papers. Therefore we have implanted 



devious traps that insure those reading this section will inevitably get sucked into reading 



the whole paper ... (just joking!). 22 
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1 Introduction. 



Reconstructing a series of 3 dimensional density distributions from a finite number of 2 
dimensional measurements is impossible unless prior assumptions of some sort are used 
fj]]. The difficulty comes from the fact that, without fairly strict assumptions, many dif- 
ferent density fields project to the same radiograph. To state this another way, a radio- 
graphic measurement device, thought of as a projection operator, has a nontrivial null 
space. Our approach in this paper is to discretize the object space and the radiograph 
(measurement) space. We then combine a sequence of radiographic measurements into one 
super- measurement. Combined with the operator which determines dynamics, the single 
time projection operator can be turned into an extended projection operator that maps 
a sequence of objects into a super-measurement. Due to the dynamical constraints, the 
dimension of the object sequence space does not grow as the length of that sequence in- 
creases. On the other hand, the size of the data space (the space of super-measurements) 
does, implying that eventually, the extended projection operator has a trivial null space. 

Now a look ahead. In section ^, we briefly outline the problem. Section || outlines the 
notation used for the rest of the paper. In section ||| we introduce the problem in it's linear 
setting. The key notion turns out to be that of transversality which we recall in section 
||. In this section, we also introduce and prove a theorem which bounds how slowly the 
dimension of the projection operator null space decreases as the number of measurements 
incorporated in the super-measurement goes up. There is a lower bound on the number of 
measurements that are needed to get a unique inversion. If q = (dimension of the object 
space) -T- (dimension of the measurement space), then the lower bound is simply \q] - the 
smallest integer greater or equal to q. We will say that a particular combination of linear 
system, L and measurement projection, P has the optimal reduction property if the number 
of measurements needed to get an unique inversion equals this lower bound. Section |6] looks 
at the prevalence of linear systems which (w.r.t. a fixed P) having the optimal reduction 
property. Next, in section [j], we outline a proof of the extension of one of the results to 
the case of nonlinear dynamics. The relation to known results is discussed in section ||. We 
look at a few numerical examples in section ^ and begin to explore the relationship between 
over-determination and noise reduction. We close with a summary and discussion. 



2 The Problem 

Radiographic experiments measuring very fast events typically produce data consisting of 
a sequence of 2-d projections. These 2-d projections are created by bombarding some 3-d 
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distribution of density - the object - with penetrating radiation of some sort such as high- 
energy x-rays or protons. The number of angles at which the data is taken is typically 1. We 
idealize this to get the following model. The object will be a point x in an object space X 
which changes from measurement to measurement as dictated by a linear operator L, acting 
on X. The measurements d which lie in the measurement space D will be generated by the 
action of a measurement or projection operator P. Thus, if the object and measurement at 
time i 6 N are denoted x% and dt respectively, we can express the actions of the operators 
L and P in the following way: xt+i = L(xt) and d t = P(xt). See figure |l[ 
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Figure 1: The Problem 

We define the extended (or experimental) spaces to be the product spaces X = X T 
and D = D T where T is the number of observation times in a particular experiment. If we 
have a particular sequence of points in the object space, then this sequence is a single point 
x, = (xi,X2, ■■■■,xt) in X. The measurement process produces a point d = (d\,d2, ...,<i*r) in 
the extended (or super-)measurement space D. This can be succinctly expressed using the 
extended projection operator P = P T since then d t = P{x t ). 

If A is defined to be the T-l by T block matrix, 
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then the null space of this operator, denoted Na, is exactly the set of elements of X 
which satisfy the dynamics. That is, Na = {x £ X\x t +i = Lx t where x = {x\,Xi, ...,xt)}- 

Let Np be denote the null space of P. The inverse problem is now solvable when 
N A n N p = {0}. 

3 Notation 

We now establish the notation we will use throughout, except in section where we find 
it more convenient to modify the notation. This section should be used as a reference. 

T = The number of radiographs. 

X = The Space of objects - we will use M. n . 

x = An element of X. 

X = X T = X x X x ... x X. 

x = An element (xi, ...,xt) of X. 

D = The space of radiographs - we will use M m . 

d = An element of D. 

D= D T = D x D x ... x D. 

d = An element (d\, dr) of D. 

L = The linear operator on X that gives the dynamics: Xj+i = L(xi). 

P = The projection (measurement) operator P : X — > D. We assume that 
P is full rank since otherwise we may choose a smaller D and consider P with this 
restricted range to get a full rank P. 

N = The null space of P. 

p = The dimension of the null space of P. 
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P = The extended or product projection operator P : X — > D. 
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A = An operator from X to X T 
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Na = null space of A 

P Na EE P(N A ). 

[A, C) = The Cartesian product A x ... x C. 

And we use the following standard notation. 

dim(i^) = The dimension of the space or subspace H 
M L = The orthogonal complement of M. 

And the almost standard notation ... 

S ffi Q = S intersects Q transversely. 

with the slight twist that when S and Q are subspaces, the intersection always includes 
the zero vector. So for example, a transverse intersection of two 1-dimensional subspaces 
of R 3 is just {0}, instead of empty intersection as would be the case for two 1-dimensional 
curves m R 3 . 
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4 Linear Stationarity implies easy solutions 



Suppose now that X = W 1 . Pick a basis for X, {bi} where i = l...n. Then Na is spanned 
by bi S X where b\ = (bi, L(pi), L T ~ 1 (bi)). For example, suppose that the linear operator 
L has a complete set of eigenvectors where i goes from 1 to n. Then Na is spanned 

by oji £ X where &i = (ui, L(oJi), L T ~ 1 (u; i )) = (uj h \u h Xj _1 u;;). 

Now suppose that D = M m . Form Er, the n by mT matrix where row i is bi = 
(bi,L(bi), L T ~ 1 (bi)). We then have that the inverse problem of getting x from P(x) has 
a unique solution iff the rank of Er is n. (This is virtually identical to the usual test for 
observability from control theory, except that I am not assuming that T = n. See for 
example || p. 178 or || p. 271). We would like to know how big T needs to be: how 
many measurements do we need? In the next section we begin to answer this question with 
an upper bound on T provided the null space of P and L do not have a special relation. 

Before we do this, let us look at Na D Np a little more carefully. If the "true" sequence 
in X is x* and we have measured d* = P(x*), then anything in N T can be added to x* 
without changing the observed data P(x*). But since we are only interested in those null 
sequences that also satisfy the dynamics, we can restrict ourselves to n* satisfying, 

n* = (m, 7*2, n T ) e[N,Nn L(N),N n L(N n L(N)),N n L(N n L(N n L(N))), ...}, 
where we are using [A, B, C] to denote A X B X ... X C. For ease of reference let 

Ni =N 

N 2 =NnL(Ni) 

N 3 =ND L(N 2 ) 

N T = N nL(N T ^). 
This permits us to rewrite the above expression for n*: 

n* = (ni,n 2 ,...,n T ) G [N u N 2 , N T }. 

Now assume that L is invertible. If Nt = {0}, it follows that E^ has rank = n and 
NA^Np = {0}. If L is not invertible, then Nt = {0} does not imply that NA^Np is trivial 
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since N A C\N p = {all possible values of(L- 1 (L- 1 (...L _1 (0)...)), -, (0), 0)}. Therefore, 
trivial iVr is a necessary but insufficient condition for the invertibility of the data d* when 
L is not invertible. 

Remark: Note that n* = (m, n 2 , n T ) G [-/V, iV n L(N), N n L(iV n L(N)), ...] is 
actually more than the set of "null" solutions. It turns out to be small enough to meet our 
needs. 

5 Transversality is (more than) enough 

We now study how the dimension of N/, depends on k. What conditions imply that even- 
tually dim(A^fc) becomes for some k? How prevalent are the L's for a fixed P that have 
T F = { minimal T such that Nt-i = {0} } = \n/d~\. That is, how prevalent are the L's 
having the optimal reduction property? 

Let us begin by reiterating the definitions of the N^. 

Ni =N 

n 2 Eivni(JVi) 

iV 3 = N n L{N 2 ) 
N T = Nr\L(N T _i). 

If G and H are linear subspaces of J then the only situation stable to small perturbations 
is that the intersection of G and H has minimal dimension. That is, we expect dim(G n H) 
= dim(G) + dim(i^) - dim(J) where a non-positive result indicate the trivial intersection 
of dimension zero. If this is the case we shall say that G and H are transverse and will 
write this as G FTT H. 

Referring to the above definitions of the iV^'s we see that these sets decrease in size at 
the maximum allowable rate if the intersections defining them are transverse. For example, 
if as above, dim(X) = n and dim(iV) = p, transverse intersections imply that the sequence 
of dimensions is p , p + p — n , p + p + p — n — n or if we note that d = n — p then the 
sequence is p , p — d , p — 2d , p — 3d and we get that dim(A^p n / d -| ) = (Remember that 
we are assuming that P is full rank.) 

Suppose that the intersections are not transverse. Then we still have the following lower 
bound on the rate at which the dimension of the N^s decrease. 
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Theorem 1 (Minimal Reduction Theorem). If there is no nontrivial subspace G of N 
such that L(G) C G then for i < dimiV + 1, dim(Ni) < dim(N)-i + 1. 

Proof. Assuming for the moment, that iVj+i C iVj, we get that dimA^ + i < dimiVj. Then, 
if dimiVj+i = dimiVj, we conclude that Ni = iVj+i. By definition of iVj+i, we then have 
N{ = N H L(Ni) so iVj C L(N{), hence dimiVj < dim L(iVj). Since L is a linear operator, 
dim L(iVj) < dimA'j holds as well, so Ni = L(iVj). But TVj is a subspace of N, so it 
follows that Ni = {0}. This means that dim(iVj + i) = dim(A r j) only if Ni = {0} so that 
dim(ATj +1 ) < dim(A r j) — 1 if dim(JVj) / 0. This implies our monotonically decreasing upper 
bound for the dimensions of the Ni. 



To see that N i+1 C A^: we use induction. We have N 2 = N n L(JV) C N = N±, so 
N 2 CNl If 7V fe+1 C AT fe , then N k+2 = N n L(A^ fc+1 ) CiVn L(iV fe ) = iV fc+1 . 



Can one find such a bad example? Yes! Consider the case in which X = E 6 and 



□ 
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p = 
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Then we get that N is given by 
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where the x's can be any value. This gives 
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' " 
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6 Optimal Dynamics are generic 

We now present and prove a theorem that addresses the point of how prevalent operators 
having the optimal reduction property are. In fact as the title of this section suggests, 
optimal L are generic. By generic we mean open and dense (not residual). 

In the following proof, we will use &a to denote dim(A). Since the set of all invertible L 
in T is open and of full measure in R( T_1 ) n we shall assume that L is invertible throughout 
the proof. We shall also use the fact that for invertible L, L(A n L{B)) = L(A) n L 2 (B). 

Theorem 2 (Extended Linear Transverse Intersection Theorem). If the set of op- 
erators L = (L\,L2, ...,Lt-i) is identified with IR( T_1 ) n and we define T C R( T_1 ) n to be 
all those L £ ]R( T_1 ) n2 such that 



Li(Ni) m N Vi. 
Then T is open and dense in R^ T ~ l ^ n ' 2 . 



Proof. We first observe that M\ (TT M 2 rfp M± (M 2 ) = min(d M ± , c?m 2 ) or equivalently that 

rank(P M ± o Pm 2 ) = min ( ii M ± ,cZm 2 )- A little bit of thought is enough to convince oneself 
that rank(PjY± o Pun)) = rank(P Ar ± o L o P/v). Therefore we get that TV HT L(N) <^ 
rank(PjY± o L o P N ) = min(d N ± , c?Ar).(Here we have used the invertibility of L to conclude 
that djy = dum). Let us approach the problem a little more generally. We shall use the 
fact that K fTT L(M) 44> rank(P ft -± o L o Pm) = min(d^±, c?m) to show that the set of 
all L in lR n such that K fTT L{M) is open and of full measure. Here K and M are linear 
subspaces of R n . 
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Define column vectors p.^ for i = 1, d K ± that are orthogonal to each other and span 
K^. Likewise let q.^ for i = 1,...,cIm be column vectors orthogonally spanning M. Define 
the n by n matrices P and Q as follows: 



and 



/ Pl,l P\,2 
P2,l P2,2 



Pl,d K± 







\ Pn,l Pn,2 ■ ■ ■ Pn,d K ± ■ • • / 



Q 



( 9i,i 91,2 

92,1 92,2 



91, d M 

92, d M 



\ 




\ 9n,l 9n,2 • • • q n ,d M • • • / 



Then 



and 



so that we get 



Now we note that 



P R± =PoP 1 



Pm = Q o Q , 



P^x oLoP M = PoP T oLoQoQ q 



rank(P o P 1 o L o Q o Q 1 ) = rank(P T o L o Q) 



To show that is open in W 1 , observe that 
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U L o 





P T o L o Q = 



v 



and therefore 



rank(P T o L o Q) = mm(d K ± , <1m) Ul is full rank. 



Note that Ul is a continuous function of U, i.e. Ul ■ IR™ 2 — > E dM ' d ^ ± is continuous 
(actually smooth). Assume without the loss of generality, that d K ± > cIm- Let (f)^ (Ul) 
be the du dimensional measure of regions in M. K± applied to the parallelepiped with edges 
equal to the columns of Ul- Then 1* is precisely equal to (C/i)^ 1 ((0 d ^ ± ) _1 (1R\ {0})). Since 

both U l and (f>^ ± are continuous, we have that is open. 

To show that R n \ % has zero n -dimensional Lebesgue measure, we first introduce a 
change of coordinates. Define P to be an orthogonal matrix obtained by filling in the zero 
columns of P appropriately. Obtain Q from Q analogously. Then 



P T oLoQ= P T oPoP T oLoQoQ T oQ 




1* ={L\ rank(L L (ul)) 



mm(d K ±,d M ) } 



={L\ L L (ul) is full rank} 
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Since P and Q are orthogonal, we have that the n 2 dimensional Lebesgue measure of 



= {L\ L L (ul) is full rank} 

and 



1* = {L\ L(ul) is full rank} 

are equal. 

We now prove that W 1 \T* has measure zero. Any L can be written as the block matrix 
with dimension of L u \ being d K ± x cIm- 



We can write L out in terms of elements as 



L 



( hi 

hi 



ll2 
h2 



\ Inl ln2 



hn \ 
hn 

Inn / 



Now assume that d K ± > d M - Identify E n2 ~ d K ^- d M+( d M^)i<i K± +i) with the e i ements Q f 
fij of an n x n matrix with the elements fd M ,d M » fd M +i,dM » ••• > fd K± ,dM removed. 
Next, define a mapping of <p : r 2 -^i"iM+(4f-i)-(«i K i+i) _^ R n 2 by 




/ij for (i, j) £ {i < d K ± and j = d M } 

E fkd M fik for (i, j) G {i < d K ± and j = d M } 



which has, as its image precisely those matrices in which L u i has column dM that is the 
linear combination of the first dM — 1 columns. If we redefine our mapping to get a series 
of completely analogous mappings, each of which has a column of L u i being dependent on 
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the other columns of L u i, then we end up with du such maps. Since each of these maps are 
smooth, and singular (the rank of the derivative is not equal to the dimension of the image 
space) Sard's theorem Q tells us that the n 2 -dimensional measure of the image of each 
map is zero. Therefore the union of the images also has measure zero. But this union is 
exactly W 1 \ T*. Since the case of d K ± < d>M is completely analogous, we have now shown 

2 

that T* is open and of full measure in M. n . 

To complete the proof we let the operator change at each step so that x 2 = L\(xi) 
, X3 = L 2 (x 2 ) j and so on. Now we have an extended operator L = [L\, L 2 , Lt-i) € 
(M n2 ) T ~ 1 . We will show that the set 6 = {L\ N^L^Ni) for i = 1,2, ...,T - 1} is open and 
dense in (R™ 2 ) 7 ^ 1 . 

2 

Define Ci to be the open subset of full measure in M. n whose members, L\, satisfy 
N ni Li(Ni). Choose a countable subset, Di, which is dense in Ci. Now for each element 
D\ of D>i, define C\ to be the open full measure subset of W 1 such that L E Cj implies that 
N m L(Nr\D$(Ni)), or equivalently, N ffi (L(N) n L ■ D$(Ni)). Define C 2 = Dit C 2- Since 

2 2 

C2 is dense in W 1 we can pick a countable D2 C C2 that is also dense in R n . We have that 
Li € Pi and L 2 G P 2 implies that N HT Li(Ni) and N (TT L 2 (N D L^N^) . Continuing this 
process we obtain Bj for i = 1, 2, T — 1 such that [L\, L2, Lr-\) G Di x I2 x ... x Dy_i 
implies that all the intersections are transverse, i.e. that iV 7TT Lj(A^) for i = 1, 2, T — 1. 
We have therefore found a subset of S which is dense in (W 1 ) T ~ 1 . 

Now we show that S is open. The requirement that each of the intersections are 
transverse is equivalent to the requirement that 

dim(iV n Li(N)) = d N + d N - n 
dim(iV n L 2 (N) n L 2 Li(N)) = M N - 2n 
dim(iV n L 3 (N) n L 3 L 2 (N) n LzL 2 Li(N)) = 4d N - 3n 

dim(iV n L T -i(N) n L T _iL T _ 2 (iV) n ... n L T - X ...L X (N)) = Td N — (T — l)n. 

which in turn is equivalent to a requirement involving orthogonal complements, specifi- 
cally that 



2(n — djv) 



riV- 1 o L^ 1 
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rN 



_L 



3(n — d 



-1 



rN 1 - o L 
rN 1 - oL^oL 



-i 

2 



T(n - d N ) 



rN 1 - 
rN 1 - o L T \ 



rN 1 o L 1 1 o ... o L T i 2 o L T i x 



all have full rank where riV^ is the matrix with rows equal to independent n-dimensional 
vectors spanning the linear subspace N ± . This last set of expressions follows from the fact 
that if M and K are linear subspaces of IR n then (M n K) 1 - = span(M ± , K 1 -) so that the 
matrices immediately above have rank 2(n — (ijv), 3(n— disf), T(n — djv) iff the the previous 
intersections have dimensions 2d^ — n, 3(ijv — 2n, ...,Tdj^ — (T — l)n respectively. But this 
last expression can be seen to be exactly those matrices which satisfy the equations 



J 2(n-d N ) 



(rows from first matrix) ^ 



^3(n-(iiv)( rows f rom second matrix) / 



^T(n-<f JV )( rows f rom T - 1st matrix) / 
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where ^(vectors) measures the j-dimensional volume of the the parallelepiped spanned 
by the vectors in M.K But since the inverse operation is continuous on the set of invertible 
matrices and these volume functions are smooth, we have that the set of (Li, L2, Lt-i) 
having the full rank property is open. Thus, & is open in (W 1 ) T ~ 1 . 

□ 

In the next section we conjecture an approach to the nonlinear case which uses the 
above derivation, but before we do this we show that the linear case can actually be made 
significantly simpler. 

Theorem 3 (Linear Transverse Intersection Theorem). If the set of operators L is 

2 2 2 

identified with W 1 and we define T C R n to be all those L € W 1 such that 

L{Ni) nr N Vi. 

2 

Then T is open and dense in R n . 

Proof. As was seen in the proof of the previous theorem, the transversality requirement and 
the fact that we are considering the case of L« = L for all i, reduces to 

n 

riV- 1 
rN 1 - o L- 1 



rN 1 - 
rN 1 - o L- 1 
rN L o L- 2 



2(n — djy) 



3(n — djv) 



T(n - d N ) 



n 

riV-L 
rN 1 - o L- 1 



rN 1 - o L- {T -^ 
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all having full rank. But this is equivalent to another full rank condition as follows. Let 
cN 1 - be the transpose of rN^. In other words, while rN 1 - are the orthogonal row vectors 
that span compliment of the null space, are the column vectors that span the "same" 
space. This condition is equivalent to the following. For a dense and open set of L: 

dim(Sfc = span(cA r± ,L o C N L , ...,L k ' 1 o cN ± )) = min(k ■ d cN ±,n) 

That is, for an open and dense set of L the sequence of subspaces generated by the 
iterates LP o cN 1 - j = 1, 2, k — 1, are of maximal dimension. 

The proof of this fact is well known. For lack of a reference I give a proof here. Without 
loss of generality let cN^ be the k left most columns of the nxn identity matrix. Then the 
dimension of Sk is the rank of the matrix formed by taking the k left columns of I, followed 
by the k left columns of L, followed by the k left columns of L 2 , and so on. 

Pick the upper left matrix minor and compute it's determinate this will give a polynomial 
in in, h2, hi, ■■■Inn which we want to show is nonzero except on an open and dense set. As 
long as the polynomial is nonzero at one point then we are done since this implies that the 
set of zeros occupies a submanifold of R n that is at most, n 2 — 1 dimensional. (So we even 
have more ... the set of "good" L has full measure and is open!) 

To show that the determinant of the upper left matrix minor is nonzero at a point, 
we consider L = permutation matrix that shifts everything to the left k clicks. This gives 
us the identity matrix in the upper left matrix minor and so the determinant in question 
evaluates to l! 

□ 

Remark: As the above proof shows, X is in fact open and full measure. This improves 
the result from one stated which implies stable approximation by L having the optimal re- 
duction property, to one that implies this AND the improbability of non-optimal reduction. 

2 

Further improvements would involve the characterization of X e = Xn {L E M. n |cond(L) < 

1/6}. 

Remark: The above proof also works with slight modifications to prove theorem ||, but 
the proof given there leads to a conjectured proof for the case of nonlinear dynamics, and 
so seems more useful even if it is more cumbersome. 
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7 Extension to the Nonlinear Case 



The above theorem is extendable to the nonlinear case as follows. Actually, we are conjec- 
turing such an extension in what follows. It should be noted that the nonlinear "extension" 
does not imply the linear theorem. 

In this section the state space (object space) will be M, a compact manifold of dimension 
m. The projection operator will be a smooth function P : M — > M. d and the dynamics will 
be given by F^s which map M to itself diffeomorphically. Instead of using linearity to get 
a fixed null space, we find the the "null" space we are now interested in is a level set of P. 
These sets can change in nontrivial ways as the point in the range of P changes. We now 
want to know about intersections of these "null" sets with images of other of the "null" sets 
under F. 

The transversality theorem found on page 74 of || implies that the set T of / € 
C r (M) which map a submanifold K of M back into M to intersect another submanifold N 
transversely, is open and dense. (For compact M the topology is nice, see [|| for details.) 
What we need is a bit more complicated. 

Let # be the quasi-stratification of M into the level-sets $ x , x £ R d of P. Let N 
be the union of a finite number of (not-necessarily injectively) immersed submanifolds of 
dimension < n whose self-intersections are transverse in the sense that the tangent spaces 
of the "participants" in the intersection span the largest possible subspace of Ti x M where 
i x is a point of self-intersection. 

Conjecture 1. For an open and dense set of F <G D°°(M) 

dim(I x ) < max(n - d, 0) Vx G M. d 
where I x = (F(N) D^ x ) and, 

I x is a finite union of stably immersed submanifolds 



This permits us to conclude that, for any initial point xq G X and dense F = (F\, 
the intersection obtained by T = \m/d~\ measurements will have dwindled to a finite set of 
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points, call it S. The next measurement (T = \m/d] + 1) will generically have precisely one 
point in the intersection of the set S and the level set corresponding to the next measure- 
ment. (We use the same argument as we used in the linear case to get a product of dense 
sets Di x ... x Dt-i-) 

Now, if indeed the \m/d\ + l"th" intersection is a single point, then the mapping 
G : x G M — > (Px, PFix, PF2F1X, PFt-i-.-Fix) has only one point in the inverse 
image of G(xq). Since the point we have "found" (the conjecture above) comes from stable 
intersections this should guarantee that the point G(xq) in fact has a neighborhood in which 
G is invertible. This should in turn guarantee that there is an open neighborhood B e of 
G in C°°(M) such that H G B t implies H*~(H(xo)) is a single point. We then use the 
fact that a small neighborhood of F maps into this small neighborhood under the mapping 
J -> {P,PJ 1 ,...,PJ T -i-Ji) for J G D°°(M). 

We have arrived at our second conjecture. 

Conjecture 2. If M is a compact smooth manifold of dimension m, P is a smooth function 
mapping M to M. d , xq is a particular point in M , T = \m/d~\ + 1, 2) = (D) T_1 is the T — l- 
fold product of the space of smooth diffeomorphisms from M to M (D°°(M) ) and we define 
H s = (P,P8i,...,P6t-i...5i) ■ M -» R dT , where 5 = (5 U ...,<5 T -i) G 2) , then the set O of 
5 G D such that 

Hf(H s (x )) = {x } 

is open and dense in D. 



8 Relation to Known Results 

The results obtained above are known as observability results in control theory and phase 
space reconstructions (delay coordinate embedings) in dynamical systems. Our results are 
different in that we consider variations of the dynamics with the observation function kept 
fixed whereas other results either assume that the dynamics are fixed and the observation 
function changes or that both the dynamics and the observation function is variable, see |p!o| , 
||, [l], ^, H . While Aeyels does consider the case where the observation function is fixed and 
the dynamics are variable he does so for vector fields (not maps). He is also looking at the 
case where he wants all initial points to be recoverable from the sequence of measurements 
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and there requires 2n+l , 1-dim measurements to recover the n-dim initial points. Similar 
comments apply to the comparison to Stark's more recent paper || . Our minimal reduction 
theorem is a more precise version of the well known theorem in control theory that states 
that if the observability matrix is not full rank then no number of measurements can give 
you full information on the state of the system and if it is full rank then you need at most 
n measurements (of any dimension) of a system that has an n-dimensional state space. 

9 Numerical Examples 

We now give two examples in which we apply the technique described in section [| to invert 
simulated sequences of (noiseless) radiographs using one view. Our purpose in this section is 
simply to demonstrate the procedure. We assume our object lies within a 10 x 10 pixelation 
and has constant density within each pixel, so the object space X is M 100 . We use the same 
initial condition with two different linear operators L\ and L 2 which we describe below. In 
each case, the projection P sums the values down the columns of the pixelation. We use, 
in a sense, the largest parameterization of our object space; namely, we assume nothing 
about the object and seek to determine the value in each pixel. At the end of this section, 
we comment briefly on the poor numerical conditioning of these problems and indicate first 
steps taken to improve the numerics. This is a subject of current study. 

To reiterate the procedure, first choose a basis {&i}J£i for X. With L and P representing 
the dynamics and projection operators respectively, we build a 100 x lOi matrix E where 
the ith row of E is (Pbi, PLbi, Pl?b{ . . . , Pl}~ x bi) = Ptbi. As soon as t is large enough so 
that rankE = 100, we have a unique solution x for the equation xE = d*. Since we know 
the dynamics L, we can then reconstruct the sequence x* = (x, Lx, L 2 x, . . . , L* x). 

In each example, we chose the canonical basis {ej}^ for X where e%{j) = 5ij (1 < j < 
100). The first linear operator L\ can be described as a combination of a diffusion and a 
shift. The effect of L\ is pictured below for various times t. 

Here, the rank of E increased by 10 each time step, so we achieved rankE = 100 in the 
minimal number of steps and were able to solve for the initial condition x. 

Our second operator, L2, was a diffusion operator where the diffusion coefficient varied 
over the pixelation. Namely, the diffusion coefficient in the i, j-pixel was |(i 3 j 2 10 _5 ) 1//4 (so 
the rate of diffusion was greatest in the lower right corner of the pixelation and was least 
in the upper left corner). The effect of L2 is pictured below for a few times t. 
In this example, the rank of E again increased by 10 each step reaching 100 after 10 steps. 
Pictured below are the initial condition x and the reconstructions obtained by using the 
data sequence (Px* , PL 2 x* , PL\x* , . . . , PL\x*) for t = 9, 14. 
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(a) X 



(b) Lfa: 



(c) Lfx 



(d) L\x 



Figure 2: Initial condition and L\x for f = 2,5,9. 
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(b) L\x 
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(c) L\x 
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(d) L\x 



Figure 3: Initial condition and L\x for t — 2, 5, 9. 

With regard to the numerical conditioning of these problems, we note that the condition 
numbers of the matrices E constructed using L\ and L2 were on the order of 10 12 and 10 11 
respectively. By running the dynamics longer than the number of time steps required to 
achieve full rank, we were able to reduce the condition number in both cases. Namely, using 
L\ for 15 time steps reduced the condition number of E to 10 11 . But with L2, using 12 time 
steps reduced the condition number of E to 10 9 , and at 15 time steps the condition number 
reduced to 10 8 . In both cases, extending beyond 15 steps gave no significant improvement. 



10 Summary and Discussion 

In this paper we examined the use of dynamics in the inversion of projection data obtained 
at a sequence of times. The main results confirm that for any fixed measurement projection 



22 



2 4 6 8 10 



2 4 6 8 10 



2 4 6 8 10 



(a) x* (b)x,t = 9 (c)x,t=U 

Figure 4: Initial condition and reconstructions using Li for t = 9, 14. 

and generic dynamics, we can simply combine the number of measurements into one large 
super-measurement which we invert to obtain the state we are trying to reconstruct. A 
following paper will deal with some aspects of the stochastic or noisy case of reconstruction 
from projections using dynamics. 

What we have established is only a first step in the direction leading to the fruitful 
combination of dynamics and measured data. Many variants of the proposed underlying 
tomography problem lead to the same abstract problem that we have begun to examine. 
For example, if the dynamical propagator f is known up to a set of parameters, the resulting 
abstract inverse problem is identical to the one stated above. There is still a state space 
one is trying to observe, only now it is n + p dimensional where p is the dimension of the 
parameter space. The projection (measurement) function now defines n+p — d dimensional 
level sets that are mapped forward by the dynamics as before. We end up needing \{n+p)/d~\ 
measurements. It seems to us that there are at least two important directions to go next. 
One is the examination of the present formulation in the presence of noise. This will bring 
us much closer to "real" situations in that the prevalence of noise makes certain problems 
which are well posed in the noiseless case, ill-posed in the presence of noise. The second 
direction is the attack of a very carefully chosen concrete problem involving dynamics that 
we understand analytically or at least numerically. This will invariably involve certain toy- 
like characteristics which should nevertheless be useful for the approach to the large, more 
realistic problems. 

Natural questions that arise include: 

• How is the problem of reconstruction from a sequence of projections related to the 
reconstruction of a 3-dim object from a spatial sequence of slices? This arises when 
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one wants to interpolate a set of CAT scans to generate a 3-dim density image. A 
related problem arises when one is trying to compress a video movie by using some 
clever interpolation in the uncompress process. 

• How do we do (algorithmically, efficiently) the reconstruction in the case of nonlinear 
dynamics. Actually, even the linear case, while conceptually simple, is not easy in the 
case of high dimension and "noiseless" (only computationally induced uncertainties 
and approximation induced uncertainties). The difficulty is that extraordinary large 
condition numbers are pervasive and so any error , like roundoff, soon overwhelms you. 
We will begin to address these issues in the next paper which looks at reconstruction 
using dynamics in the presence of uncertainty. 

• If one has a set of measurements, how does one use the knowledge of the underlying 
dynamics (incomplete maybe) and the freedom to choose the object (state space) 
parameterization in such a way as to get a well posed inverse problem with as little 
as possible "wasted" measured information. That is, How does one use all the prior 
information and measured information in the generation of the final reconstruction. 
Even if we are using all our information, there are different ways of distributing 
remaining uncertainty about the reconstructed object. At each for a given data 
set what sort of different parameterizations give rise to this no wasted information 
situation? 

• In the high dimensional case, the questions asked in the previous bullet are incredibly 
hard to answer. What approximate answers can be generated? Can we tell how far 
from the optimum we are? for example, can we obtain bounds on the amount of 
information that our parameterization/reconstruction/use-of-priors wastes? 

• Suppose we do the whole analysis with e fattened null spaces. In this case, what sort 
of volume do we get for the final intersection (which before was just one point)? This 
is along the same lines as the first remark at the end of section ||. 
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